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The probability current statistics of two-dimensional open 
chaotic ballistic billiards is studied both analytically and nu- 
merically. Assuming that the real and imaginary parts of the 
scattering wave function are both random Gaussian fields, we 
find a universal distribution function for the probability cur- 
rent. In by-passing we recover previous analytic forms for 
wave function statistics. The expressions bridge the entire re- 
gion from GOE to GUE type statistics. Our analytic expres- 
sions are verified numerically by explicit quantum-mechanical 
calculations of transport through a Bunimovich billiard. 
PACS 05.45.Mt,05.60.Gg,73.23.Ad 



I. INTRODUCTION 

For quantum chaotic closed system it is well known 
that the statistical properties of the energy levels are 
described by random matrix theory (RMT) Q. They 
follow the Gaussian orthogonal ensemble (GOE) and 
the Gaussian unitary ensemble (GUE), depending on 
whether time-reversal symmetry (TRS) of a system is 
preserved or not. In the same way the wave function 
statistics obeys different laws in the two cases. Let the 
scaled local density be p(r) = A\ip(r)\ 2 where ^(r) is 
the normalized wave function and A is the area (vol- 
ume) of the system. As prescribed by GOE the proba- 
bility distribution is the well-known Porter-Thomas (P- 
T) distribution P(p) = (l/y/Wp) exp(-p/2) when TRS 
is present (the Hamiltonian H is invariant under time 
reversal T(t — ► — t) and ip may be chosen real). On 
the other hand, the distribution takes the exponential 
Rayleigh form P(p) = exp(— p) as described by GUE 
when TRS is broken (H is not invariant on T and yb must 
be complex). It is easy to understand qualitatively why 
the statistics is so different in the two cases, for example, 
why small values of p have a much larger weight in GOE 
than in GUE. In the first case the real wave function van- 
ishes along nodal lines in two-dimensional (2D) systems 
(surfaces in 3D). On the other hand, yb is complex in the 
second case and vanishes only at nodal points (lines in 
3D) resulting in less probability for small p. Depending 
on the relative weights of the orthogonal real and imagi- 
nary parts of i/j one can also define intermediate statistics 
that applies to the entire crossover region from GOE to 
GUE || . 

We now consider what happens when the system is 



made open, for example by attaching electron leads to 
some exterior reservoirs and a stationary current through 
the system is induced by applying suitable voltages to 
the reservoirs. The additional flexibility gained in this 
way leads to a number of interesting cases for the wave 
function statistics. Let us first look at the case when 
there is no current flow. The statistics will then be the 
same as for the closed system above, i.e., the kind of 
statistics simply depends on whether the Hamiltonian is 
invariant under T or not. On the other hand, if there is 
a stationary current via the leads we have to deal with a 
scattering wave function. This function, which must be 
complex, is written in 2D as 



yj{x,y) = u(x,y) + iv{x,y). 



(1) 



and satisfies (V 2 + k 2 )yb = 0. Even if the Hamiltonian 
itself is invariant under T, the statistics will not follow 
GOE since the scattering wave function is complex. Be- 
cause of the boundary conditions associated with the 
scattering wave function it is no longer an eigenstate of 
the usual time reversal operator CT where C is complex 
conjugation of ip(r). We will show how the two Gaus- 
sian random fields u and v are identified to recover the 
intermediate statistics discussed above and how a univer- 
sal distribution for the probability current density can be 
found. We will also compare theory with explicit numer- 
ical calculations for an open Bunimovich stadium. 



II. THEORY 

In the following derivation of the wave function and 
probability current statistics we assume that the real and 
imaginary parts of yb can be viewed as two independent 
isotropic Gaussian fields. An explict example of such 
a state is given in ref. p| in the form of a Berry-type 
wave-chaotic function. In general the assumption of in- 
dependent fields can only make sense if we first extract 
a common phase factor. This feature will turn out to be 
most useful. Let us introduce the notations 
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(uv) = 7, 



<H = < + < = <M >> (u)=0, (v) = 0. 
We define the averages as 
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where A is the area to be sampled. In our case it will 
be area of the cavity, but in principle it could be any 
area that one may wish to diagnose. In what follows we 
assume that wave function ip(r) is normalized as 



f d 2 r|V(r)| 2 = l, 

J A 



(5) 



and therefore a 2 A=l. To bring ip to a "diagonal" form 
in which the real and imaginary parts are independent 
Gaussian fields we introduce the new functions p(x, y) 
and q(x, y) by changing the phase as 



ip(x, y) -> e m ip(x, y) = P{x, y) + iq(x, y). 



The condition (pq) = now allows us to determine a. 
By this step we will also be able to find the analytic 
expressions for the wave function and probability current 
statistics. Straightforward algebra gives 



tan 2a = 



2l 



1 



a 2 - y/o*-^olal - 7 2 ) 



(7) 



Next let us consider the cumulative distribution G(p) for 
the scaled density p(r) — A\ip(r)\ 2 : 



G( P ) 



c(p) 



f(p, q)dpdq. 



(8) 



The integration is defined by the circle C (p) in the (p, q)- 
plane centered at the origin and with radius yp(r), i.e., 
(p 2 + q 2 )/cr 2 < p(r) in the integral above. The function 
f(j>, q) is the joint distribution for the random Gaussian 
fields p and q: 



f(p,q) = 



After integration of Eq. (||) we obtain 



(9) 
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1 — exp[— p(J,(n + vcosO)] 

fJL + V COS 6 



d9, (10) 



where we introduced the following notations 
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Differentiating (|l(]) with respect to p we find the final 
expression for the density distribution 



(12) 



where Io(z) is the modified Bessel function of zeroth or- 
der. 

The distribution in Eq. ( |l2| ) coincides with the results 
obtained from RMT for closed systems |§,| and is there- 
fore not new. For weakly open systems with point con- 
tacts Seba et al. || have related the statistical properties 
of the scattering matrix elements with the distribution 
P{p) and have obtained the expression above. Analo- 
gous derivation of ours is also found in Ref. [Q. Our 
way of deriving Eq. (|l^), however, explicitly shows how 
to identify the two independent random fields in a given 
wave function. For the random wave function in Eq. 
(^) the limits {q 2 ) — » 0, e — > correspond to a closed 
system and as a consequence one recovers the P-T dis- 
tribution and GOE statistics. On the other hand, the 



(6) case (p 2 ) -> (q 2 ), 



1 corresponds in this context to 



an open system through which there is a current flow. 
Consequently one finds the exponential Rayleigh distri- 
bution that corresponds to the GUE statistics. In the 
crossover region, the value of e is obtained numerically 
using Eqs. (g), (||), (0), and (jll|), i.e., e is not merely 
a fitting parameter. We have recently verified this type 
of crossover in wave function statistics for a Bunimovich 
stadium using numerical scattering methods ||. 

In view of all previous work on the generic form of wave 
function statistics in chaotic systems it is surprising that 
no attention has been paid to the corresponding current 
distributions except Ref. where they only relate the 
average of squared current to (p). Since currents may be 
measured it is of interest to establish a form also for 
currents. Below we will show how to find a useful form 
that is both simple and universal. Let us limit ourselves 
to the case of a weak net current between narrow input 
and output leads. Inside the cavity, however, there will 
be a rich, whirling flow pattern, which is strongly influ- 
enced by the vortical motions around the nodal points 
associated with the complex form of the wave function. 
Hence the net current through the through the billiard 
turns out to be only a tiny fraction of the total inter- 
nal flow, in particular for asymmetric arrangements of 
leads and wave lengths that are small compared to the 
dimensions of the cavity. As a result the corresponding 
distributions may to a good approximation be chosen to 
be isotropic. Hence the components of the current ef- 
fectively average to zero. These assumptions arc verified 
by the numerical calculations to be discussed in the next 
section. 

Our complex wave function (Q) carries the probability 
current density (h = m = 1) 



j = Im(?/>* Vi/") = pS/q — qS/p. 



(13) 



To find the corresponding distribution it is convenient 
to begin with a characteristic function for the compo- 
nents of the probability current density 



0(a) = (e la ' j ) = (exp[i(pa -Vq-qa- Vp)]>. (14) 

Since (p\7q) = (q\7p) = for isotropic fields Vp and 
Vq are statistically independent of p and q. They have 
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the same distribution as in Eq. (Q) with dispersions 
((Vp) 2 ) = k 2 (p 2 ) and ((V<?) 2 ) = k 2 (q 2 ) which fol- 
lows from the Schrodinger equation. Using the relation 
((aVp) 2 ) = a 2 k 2 (p 2 )/2 and similarly for Vq we obtain 



9(a) 
where a = lal and 



1 



1 + t 2 a 2 



= k 2 (p 2 )( q 2 )/2. 



(15) 



(16) 



From Eq. (|T^) it is now easy to calculate the distribution 
functions. For the components we have 



2n 



e{\a x \)e- la ^ da x 



1 (17) 
2^exp(-\jx\M 



and the same for P(j y ). In order to derive the distri- 
bution function for the absolute value of the probability 
current density let us consider the joint distribution func- 
tion 



p (jxJ y ) = ^-l a,J (aj)e(a)da 



1 
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(18) 



where j = |j| and K (z) is the modified Bessel function 
of the second kind. Since this expression is radially sym- 
metric one can find the probability density function P(j) 
for j by just multiplying Eq. (18) with a factor of 2nj. 
This gives us the final expression 



(19) 



III. NUMERICAL RESULTS 

As a numerical verification of the analytic expressions 
for the probability current distribution, we consider an 
open 2D Bunimovich hard-wall stadium (see inset in 
Fig.[TJ). It is characterized by the radius of a semicircle 
a and the half-length of a straight section I, and cou- 
pled to a pair of leads with a common width d. Here we 
choose (a = I) and (d/y/A = 0.0935) for which the bil- 
liard is maximally chaotic and weakly open, respectively. 
To find the scattering wave function for particles entering 
and leaving the cavity via the leads we solve the time- 
independent Schrodinger equation for ip under Dirich- 
let boundary conditions using a plane-wave-expansion 
method |?]], which gives reflection and transmission am- 
plitudes for a given energy. The wave functions are used 
to compute the different parameters entering the statis- 
tics using the explicit expressions stated above. 

Figure |l| shows the transmission probability T as a 
function of Fermi wave number k for an incoming wave 



with transverse mode n in the leads. There is a sequence 
of overlapping resonances which become broader in the 
high energy region shown in the lower section of Fig. [j]. 

For the statistical analysis of the scattering wave func- 
tions we select two typical cases: (A) A low energy with 
only one fully open channel (n — 1) for which T reaches 
unity; (B) A high energy with n = 4 and an intermedi- 
ate value for T. For the statistics the spatial average is 
taken over the billiard region corresponding to the closed 
stadium. For convenience this area is set equal to unity. 
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FIG. 1. Transmission probability T as a function of Fermi 
wave number k for the open stadium billiard: (A) A low en- 
ergy case with one open channel in the leads (n — 1), (B) A 
high energy case with n = 4. The inset shows the hard- wall 
Bunimovich stadium and the positions of the leads. 

Figure |2| shows the numerical results for P(j), P(jx) 
and P(jv) together with the analytical predictions in 
Eqs. ( |l^ ) and (|l9|). In the case A there is almost no 
reflection and hence the system is completely coupled to 
the open channel. The current statistics shows, how- 
ever, that e = 0.32, i.e., intermediate between a closed 
and fully open cases. Nevertheless, the numerical results 
show good agreement with the theory. 

Also in the high energy region B in Fig. the probabil- 
ity current distributions are well described by the theory 
as shown in Fig. ||[ Here e = 0.86 which is close to the 
exponential Rayleigh (GUE) case e = 1. The transmis- 
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sion is, however, lower than in the previous example. As also the case of thin microwave resonators |Ip,|iO|] and 
it appears from these there is no simple relation between applying the present theory to the Poynting vector. 
T and e. 
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FIG. 2. Distribution of probability current density P(j) 
(top) and its components P(jx) and P(j y ) (middle and bot- 
tom) in the open stadium billiard for the case A in Fig. |l| 
Solid curves show the analytical predictions for e = 0.32. (For 
convenience K = l.m = 1.1 



IV. CONCLUDING REMARKS 

We have derived the statistical distributions for wave 
functions, probability current densities and correspond- 
ing components for 2D open quantum systems with clas- 
sically chaotic dynamics. The expressions for the proba- 
bility currents are universal in the sense that the shape of 
the distributions is independent of the mixing parameter 
e, i.e., only the width changes with e. This is in contrast 
to the wave function statistics that transforms gradually 
from GOE to GUE type with increasing e. Obviously 
these ideas carry over into 3D. 

The statistics for mesoscopic transport through a 
chaotic open billiard was also studied numerically with 
enough statistical resolution to compare with the analyt- 
ical predictions. The discussed results give a numerical 
proof of the predictions for the probability current den- 
sity developed here. It also appears that experimental 
verifications are possible. For example, images of the 
coherent electron flow through a quantum point contact 
have been observed in recent experiments H. There is 
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FIG. 3. Same as in Fig. |^ but for the case B for which 
e = 0.86. 
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